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We study a Bose-Einstein condensate (BEC) of a dilute gas with dipolar interactions, at finite 
temperature, using the Hartree-Fock-Bogoliubov (HFB) theory within the Popov approximation. 
An additional approximation involving the dipolar exchange interaction is made to facilitate the 
computation. We calculate the temperature dependence of the condensate fraction of a condensate 
confined in a cylindrically symmetric harmonic trap. We show that the bi-concave shaped con- 
densates found in Ref. [lj in certain pancake traps at zero temperature, are also stable at finite 
temperature. Surprisingly, the dip in the central density of these structured condensates is actually 
enhanced at low finite temperatures. We explain this effect. 
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I. INTRODUCTION 

The realization of a Bose-Einstein condensate (BEC) 
of 52 Cr Q marked a major development in degener- 
ate quantum gases in that the inter-particle interac- 
tion via magnetic dipoles in this BEC is much larger 
than that in alkali atoms, and leads to an observ- 
able change in the shape of the condensate. The long 
range nature and anisotropy of the dipolar interaction 
pose challenging questions about the stability of the 
BEC and have led to predictions of unique phenom- 
ena, such as roton-maxon spectrum, different phases 
of vortex lattices, and bi-concave shaped condensates 
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The effect of finite temperature on these phenomena, 
or new temperature-dependent effects, remain largely un- 
explored. Theoretical finite temperature studies have 
been confined to path integral Monte Carlo simulations 
Q of a small number ( f 00) of particles, and to a ho- 
mogeneous, quasi ID system (l3| . In the latter work, 
the Popov approximation to the HFB theory has been 
applied [3, fl5| . In systems with short range interac- 
tion, the Popov approximation has been found to give 
excitation spectra in good agreement with experiment 
for temperatures up to half the critical temperature for 
condensation [r|. For the density profile, good accu- 
racy was shown for even higher temperatures, up to the 
critical temperature. 

In [l2| we have introduced a new computational al- 
gorithm which allowed us to calculate the Bogoliubov 
excitation spectrum of dipolar condensates in cylindri- 
cally symmetric 3D traps, at zero temperature. Here 
we extend this work in a natural way to describe finite 
temperature properties in the Popov approximation. In 
essence, the Popov approximation is a self-consistent so- 
lution in which the Bogoliubov excitation spectrum is 
computed and the different excitation modes populated 



according to Bose statistics. This leads to depletion of 
the condensate, and thus a shift in excitation frequencies. 
The excitation spectrum is re-calculated iteratively until 
self-consistency is achieved. 



II. FORMALISM 

The HFB-Popov equations for the case of short range 
interactions have been described by Griffin 14]. The 
generalization to long range interactions is straightfor- 
ward, and we therefore briefly formulate it in this section. 
The confined Bose gas is portrayed as a thermodynamic 
equilibrium system under the grand canonical ensemble 
whose thermodynamic variables are the temperature T 
and the chemical potential \x. There is a one-to-one rela- 
tionship between the chemical potential \x and the total 
number of particles N. Below the critical temperature, 
Nq of the atoms are in a condensate state. The system 
Hamiltonian K for the system with a fixed chemical po- 
tential is then obtained from the Hamiltonian H of the 
system with fixed number of total particles via a Legen- 
dre transform, and has the form 

K = H — fiN = J dr&(r)(H - fJt)V(r) 

+ 1 J &( r )&(r')V(r' - r)*(r')*(r), (1) 

where $(r) is the Bose field operator that annihilates 
an atom at position r; iJ = (—h/2M)V 2 + V tra p(r) 
contains the kinetic energy and the trap potential, and 
V(r) is the particle-particle interaction potential. For 
the system treated here, the trap potential is cylindri- 
cally symmetric: Vt ra p(r) = M(u 2 p 2 +w 2 z 2 )/2, where 
M is the atomic mass, and to p and u> z are the radial and 
axial trap frequencies, respectively. For polar gases the 
potential V(r) may be written: 
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where a is the scattering length, d the dipole moment, 
r the distance between the dipoles, and 9 is the angle 
between the vector r and the direction of polarization, 
which is aligned along the trap z-axis . 

The Bose field operator is decomposed into a c-number 
condensate wave function plus an operator describing the 
non-condensate part, *fr(r) = \fTT^(j)(r) + ip(r), and in- 
serted into Eq. |T]). Terms cubic and quartic in ip{ r ) 
are treated within the mean-field approximation and the 
grand-canonical Hamiltonian reduces to a sum of three 
terms: K = Kq + K\ + the first term Kq is a c 
number, the second term, K\, is linear in tp(r) (and its 
hermitian conjugate), and the last term, if 2, is quadratic 
in these quantities. Within the Popov approximation, 
the so called anomalous terms arising from mean field 
averages of the form < tpip > are ignored and only 'nor- 
mal' terms of the form < 7/^-0 > are included [l4j. For a 
system in equilibrium the linear term K\ is required to 
vanish identically, giving a generalized Gross-Pitaevskii 
(GP) equation for </>(r): 



H + / dr'{N Q \<p{r')\ 2 + n{r')V{r' - r) 



(r) 



+ / »(r', r)V{r' - r)<j>(r') = ^(r), (3) 



where h(r) = (ip'(r)tp(r)) is the density of the non- 
condensate (thermal) atoms, and n(r' , r) = (tp'(r')^(r)) 
is the 1-particle reduced density matrix, or the correla- 
tion function, of the non-condensate atoms. The term in- 
volving h(r) represents the mean field contribution due to 
direct interaction between thermal cloud and the conden- 
sate. The term involving n(r' , r) represents the contribu- 
tion of exchange interaction between the thermal cloud 
and the condensate. Note that for V with only short 
range interaction, the exchange term reduces in form to 
that of the direct one. But for long range interaction, 
the non-local correlation function is needed. Also, note 
that a more careful treatment [UGll reveals that, for an 
atom number conserving system, the factor Nq in Eq. ([3]) 
should be replaced by Nq — 1 . This correction is negligible 
for N >> l. 

The term K 2 has the form 

K 2 = j drft (r)£ij)(r) 

+ / drdr'ft(r')n(r, r')V(r' - r)^(r) 



+ — drdr'i)\r')ft{r)V{r' - r)<p(r)<f>(r') 

+ ~T J drdr '^( r ')^( r ) v ( r ' - r ) ( t ) *( r ) ( t ) *( r ')( 4 ) 

where £ = Ho — fi+ J ' dr'V(r' — r)n(r'), in which n(r') = 
n(r') + Nf)\c/)(r')\ 2 is the total density, and n(r,r') = 
h(r, r') + Nq4>* (r)(f>(r') is the total correlation function. 



The term K2 can be diagonalized by the Bogoliubov 
transformation, 



(5) 



if the quasi-particle amplitudes Uj(r) and Vj{r) satisfy 
the coupled HFB-Popov equations: 



EjUj(r) = Cuj(r) + J dr'V(r' — r)n(r' ,r)uj(r') 

+ N [ dr'cf>(r')V{r' - r)v{r')<j){r), (6a) 



E iVi {r) = Cv j (r) + Jdr>V { r>-r)n { r,r>)v j (r') 

+ N dr' <t){r')*V{r' ~ r)u(r')4>{r)* . (6b) 



The quasi-particle annihilation and creation operators 
ctj and ctj satisfy the usual Bose commutation relations. 

In terms of it's and f's, the thermal density correlation 
function is written as: 

n(r',r) - ]T [u* (r') Uj (r) + Vj (t>J (r)] N^E,) 
j 

+v j (r')v*(r), (7) 

where N ex (Ej) is the Bose distribution for the quasi- 
particle excitations: 



N ex {E 3 



I 



exp 



\k B T 



- 1 



(8) 



The expression for n(r) in terms of u, v is obtained by 
setting r' = r in Eq. 0. Similar expressions may be 
easily obtained for the total and release energy of the 
condensate in the trap. 

As usual, the self-consistent HFB-Popov equations © 
and ([6]) are solved iteratively. At first, one sets n(r', r) = 
n(r) = 0. The thermal component contribution is then 
updated at each further step using Eq. ([7]) until conver- 
gence is reached. 

It can be appreciated that long range interactions 
present a significant challenge to the computational im- 
plementation of the HFB-Popov method. The difficul- 
ties arise even for zero temperature, where (ignoring the 
negligible quantum depletion) the HFB-Popov equations 
reduce to the Bogoliubov-De Gennes (BdG) equations, 
which due to the long range exchange interaction, are 
now integro-differential rather than simply differential 
equations. In Ref. [l2j we have introduced a new algo- 
rithm, which enabled us to solve the BdG equations for 
a gas with dipolar interactions in a 3D trap with cylin- 
drical symmetry, by utilizing the cylindrical symmetry to 
reduce the effective dimensionality of the problem from 
3D to 2D. For the finite temperatures with which the 
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HFB-Popov method is concerned, the long range inter- 
actions introduce an additional difficulty: Eqs. © and 
(J6j) involve not only the thermal density h(r), but also the 
thermal correlation function n(r',r). The most difficult 
terms are those involving the long range exchange inter- 
action such as J dr'V(r' — r)n(r',r)uj(r') in Eq. ([6]). In 
the case where there is no thermal component, the to- 
tal correlation function is expressed as a direct product: 
n(r',r) — <fi* (r')(f>(r). In this case, the above exchange 
term may be effectively evaluated by the use of Hankel- 
Fourier transform [l2| . In the presence of a thermal com- 
ponent, the total correlation function does not have such 
a direct product decomposition. It might be possible 
to circumvent this complication by direct evaluation of 
the exchange integral on a spatial grid (without resort to 
Fourier transforms; see for example [19() - however, this 
is complicated by the singular nature (1/r 3 behavior) of 
dipolar interactions at the origin. 

To obtain a feasible numerically solvable problem, we 
therefore make an additional approximation: in Eqs. ([3]) 
and (|6|) we let n(r',r) — for r' ^ r, or equivalently, 
we let n(r',r) = cf>(r')*4>{r). Physically, in Eq. ^ this 
amounts to ignoring the forces on the condensed part due 
to the long-range exchange interaction with the thermal- 
ized part. In Eq. ^ this amounts to ignoring the long 
range exchange interaction between the thermal compo- 
nent of the gas and itself. On the other hand, in Eq. ^ 
we do take into account the effect on excitation modes 
and frequencies due to the the long range exchange in- 
teraction with the condensate. We treat exactly the long 
range direct interactions which involve n(r). We also 
treat exactly the short range (contact) interaction, for 
which the exchange and direct terms are identical. 

A partial justification for the above scheme may be 
found in the good agreement between the "two-gas 
model" of dilute BECs and the full HFB-Popov descrip- 
tion for gases with short range interactions [20j. In the 
two gas description, the condensate wavefunction is that 
of a T — BEC with the appropriate (depleted) num- 
ber of condensate atoms No, and the surrounding ther- 
mal cloud is described by the statistical mechanics of 
an ideal gas in the combined potentials of the trap and 
the cloud-condensate interaction. The full HFB-Popov 
description may be reduced to the two gas description 
by letting n(r, r') = for all r, r', so that, in particular, 
n(r) = n(r, r) = 0. The reason for the success of the two- 
gas model seems to derive from the fact that the thermal 
component is typically much more dilute than the con- 
densate part. Therefore, to a good approximation, it 
may be described as an ideal gas. The approach sug- 
gested above for the treatment of dipolar BEC at T > 
may be described as treating the thermal component as 
'partly ideal', i.e, 'ideal' only with respect to long range 
exchange interactions. Thus, it is a compromise between 
the full HFB-Popov method and the two gas description. 

Moreover, We note as a general thermodynamical 
property that the correlation function n(r',r) naturally 
decreases towards zero with increase in temperature (for 



r' 7^ r). Thus, with increasing temperature, it makes 
sense to ignore the thermal long range exchange inter- 
action which is due to the correlation function of the 
thermal component of the gas. 

The number of excitation modes that need to be taken 
into account in Eq. © in order to saturate the ther- 
mal cloud density is very large (tens of thousands). For 
the higher excitation modes, the semi-classical descrip- 
tion has proved very useful and accurate [H, [2l[. Thus, 
we follow the approach of solving Eqs. (O for the dis- 
crete modes up to an appropriate energy cutoff, and us- 
ing their semi-classical version [12| for modes above this 
energy cutoff. The energy cutoff is typically somewhat 
larger than the chemical potential and is adjusted in each 
specific case until convergence is achieved. 

III. RESULTS 
A. Cr in a pancake trap 

We first study the effects of temperature for a 52 Cr 
gas in a trap with frequencies lo p = 2tt x 100Hz and 
u> z = 2ir x 400iJz (a pancake trap). The magnetic dipole 
moment of polarized Cr is relatively large for atoms, 6 
Bohr magnetons. However, the resulting dipole-dipole 
interaction is still small compared to the strength of the 
short range interaction (scattering length a = 96ao |22|). 
A useful parameter here is the dimensionless quantity 

,2 

( dd = ^2/ a ■ A homogeneous condensate is unstable if 
£dd > 1 [23] • For Cr, tdd = 0.16. However, using a Fes- 
hbach resonance, it is possible to reduce the scattering 
length and thus increase the dipolar effects (T. Pfau, pri- 
vate communication). For the present study we assume 
a reduced scattering length of 20a , so that e^d = 0.8. 
The number of atoms in the trap is taken to be 10 . 

In Fig. ([1]) we show the condensate fraction in the pan- 
cake trap as a function of temperature. For comparison, 
we have included results of an un-polarized gas by set- 
ting the dipole moment d — 0. It is seen that the effect 
of the polarization of the gas is to decrease the conden- 
sate fraction at any given temperature, compared to the 
non-polarized gas. As a result the critical temperature 
is also reduced. This effect may be expected since in a 
pancake trap the average dipolar interaction is repulsive, 
thus the thermal effect due to polarization is similar to 
that of increasing the scattering length. 

In Fig. (fj)) we plot the eigenfrequencies of the lowest 
collective modes. Although the dipolar interaction causes 
large shifts in the frequencies, the temperature depen- 
dence of these shifts is very small, except very near the 
critical temperature. The shift in frequencies due to po- 
larization is, qualitatively, similar to that of increasing 
the effective short range repulsion. Note the Kohn mode 
(uj/uj z — 1) which should remain constant at 1. The 
slight deviation from 1 at higher temperatures is due to 
Popov approximation which in effect computes the dy- 
namics of the condensate in the presence of a static ther- 
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FIG. 1: Condensate fraction as a function of temperature 
in a pancake trap with aspect ratio 1:4 (i.e, uj p /uj z = 1/4 
): un-polarized, i.e, dipole moment set to zero (dashed line), 
polarized (solid line), and ideal gas (dotted line). The total 
number of 52 Cr atoms is 100000. 
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FIG. 2: Lowest excitation frequencies a function of tempera- 
ture for a 5 2Cr condensate in a pancake trap with aspect ratio 
1:4; polarized (solid lines)and un-polarized (dashed lines). For 
each angular momentum number m = 0, 1, 2, 3, we plot only 
the lowest mode with this m. 

mal component. A more correct description should treat 
both components dynamically. Nevertheless, the devia- 
tion of the Kohn mode from the theoretical value of 1 is 
small. 



B. Cr in a cigar trap 

Figs. ([3]) and ([4]) show the condensate fraction and low- 
est collective mode frequencies as a function of tempera- 
ture for a 52 Cr gas in a cigar trap with uj z = 2ir x 100 Hz 
and ojp — 2ir x 400 Hz (with a reduced scattering length 
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FIG. 3: Condensate fraction as a function of temperature in 
a cigar trap with aspect ratio 4:1; un-polarized (dashed line) 
polarized (solid line), and ideal gas (dotted line). 
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FIG. 4: Lowest excitation frequencies as a function of temper- 
ature for a 5 2Cr condensate in a cigar trap with aspect ratio 
4:1; polarized (solid lines) and un=polarized (dashed lines). 
Only modes which are even in the z direction are shown. For 
each angular momentum number m = 0, 1, 2, 3, we plot only 
the lowest even mode with this m. 



a = 20ao as before). Notice that now the effect of the 
dipolar interaction is to increase the condensate fraction 
at any given temperature, thus increasing also the critical 
temperature for the onset of condensation. Again, this 
can be understood due to the dipolar interaction being 
effectively attractive in a cigar geometry. Similar to the 
the case of a pancake trap, the dipolar interaction leads 
to significant shift in the frequencies of the low modes, 
but these shifts depend only weakly on temperature. 
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FIG. 5: Condensate fraction as a function of temperature in a 
pancake trap with aspect ratio 1:7, for which bi-concave struc- 
ture is formed at T=0; un-polarized (dashed line) polarized 
(solid line), and ideal gas (dotted line). The total number of 
52 Cr atoms is 16300. The scattering length is 0. 



C. Bi-Concave condensates 

We now turn to examine the finite temperate effects on 
the bi-concave shaped condensate reported in Ref. 
There, we found an interesting novel structure of pure 
dipolar condensates in pancake traps at zero tempera- 
ture. For appropriate choice of parameters, the conden- 
sate density does not obtain its maximum in the center of 
the trap. Rather, the maximum density is obtained along 
a ring, and the center of the trap is local minimum of the 
density. This gives rise to a bi-concave condensate shape 
similar to that of a red-blood cell. Recently, other shapes 
have been predicted in non-cylindrically symmetric traps 
[24| . In this section, we investigate the temperature effect 
on the bi-concave condensate in a cylindrically symmetric 
trap. 

In Fig. ([5]) we plot the condensate fraction as a function 
of temperature for a pure 52 Cr dipolar condensate (i.e, 
where the scattering length has been tuned to zero via a 
Feshbach resonance), in a pancake trap with aspect ratio 
<jjp/<jj z — 1/7. For number of particles N — 16300, a 
bi-concave structure is formed at T = 0. In this figure 
it is notable that the dipolar interaction brings about 
significant change (about 10%) in the condensate fraction 
for tempeartures of order half the critical temperature. 
Yet, the critical temperature itself is almost un-changed. 
Indeed, the analytical formula of Refs. [25|,[26| predicts a 
very small reduction of 0.5% in the critical temperature. 

An interesting question is: what happens to the bi- 
concave structure with the increase in temperature? The 
Bi-concave shape eventually washes out by T — T c , yet 
(as we shall show) the shape persists to surprisingly high 
temperature. To study this, We define a contrast param- 
eter, c = 1 — n(0)/n max , where n(0) is the central density 
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FIG. 6: Biconcave contrast as a function of temperature in 
a pancake trap with aspect ratio 1:7. Solid line: contrast 
parameter for the total density; Dashed line: contrast param- 
eter for the condensed part. Inset: illustration of a typical 
biconcave density profile showing how the biconcave contrast 
c is defined. 



and n max is the maximal density. For a normal density 
profile where the maximal density is obtained at the cen- 
ter, c = 0. In Fig. ([6]) we plot the biconcave structure 
parameter for the total density profile, as well as for the 
condensed part alone, as a function of temperature. It is 
seen that when the temperature approaches about 70% 
of the critical temperature, the biconcave structure dis- 
appears. Generally, one would expect the disappearance 
of the biconcave structure due to the thermal excitations. 
We note, that according to Ref. the biconcave con- 
trast (at T = 0) is reduced with decreased number of 
particles. Thus, when the condensate is depleted, we 
also expect the bi-concave parameter to decrease. For 
T > 0.15T C , this is indeed the case. But for lower tem- 
peratures, we see that the biconcave contrast, for both 
the total density and the condensed part alone, slightly 
increases with temperature. 

To understand this effect, let us first consider the den- 
sity profile of the thermal cloud alone. Consider first 
the simplest case of an ideal gas in a harmonic trap. The 
thermal cloud occupies harmonic oscillator states accord- 
ing to Bose statistics. The lowest and most populated 
excited state, one above the ground state, has a node at 
the center of the trap. Thus, at low temperatures, one 
expects the thermal cloud to have reduced density at the 
center of the trap, even in the absence of repulsive short 
range interactions. This effect is easily verified by numer- 
ical simulations, and is seen to be more pronounced when 
the dimensionality is reduced, such as in highly pancake 
or cigar traps. Of course, for an ideal gas, the total den- 
sity (thermal+condensate) still has its maximum at the 
center. 

Consider now the dipolar gas. As the temperature is 
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FIG. 7: Radial density profiles of the condensate density (solid 
line) and lOx the thermal component density (dashed line). 
Note that we scaled the thermal component density by a fac- 
tor of 10 for visual comparison. The trap aspect ratio is 1:7 
as in Fig. [6] and the temperature is T = 0.2T C , where T c is 
the critical temperature. 



FIG. 9: Lowest excitation frequencies as a function of tem- 
perature for a 5 2Cr condensate in a pancake trap with aspect 
ratio 1:7. Solid line: lowest m = mode;dashed line: m = 1 
mode; dash dotted line: m = 2 mode; dots: m = 3 mode; 
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FIG. 8: The mean dipolar field due to the thermal component 
density whose radial profile is shown at Fig. [7] 



raised from T = 0, the lowest and most populated ther- 
mal gas mode is a sloshing (Kohn mode) with a node 
in the center. Thus for low temperatures the thermal 
component has lower density at the center, as is demon- 
strated in Fig. ([7]). The thermal component creates in 
turn a mean field, which is shown in Fig. (|8|). It shows 
that the maximum mean field potential is obtained in 
the center of the trap, even though this is not where the 
place of maximum thermal component density. This is 
due to the long range nature of the dipolar interaction: 
the contributions from the ring of maximal density of the 
thermal cloud add together in the center of trap. The 



mean field due to the thermal cloud causes the conden- 
sate part to be repelled from the center. The result is that 
the biconcave contrast of the condensate is larger than 
it would have been with the same number of condensed 
particles in a pure harmonic trap with no thermal com- 
ponent. This gives rise to the behavior seen in Fig. ^ 
at low temperatures. 

A caveat is that, as mentioned above, we made an ap- 
proximation in our computations by ignoring thermal ex- 
change effects. This may, in principle modify the effect 
of the thermal component on the condensed part, at low 
temperatures in particular. Exploring the effect of the 
thermal exchange interaction would require a consider- 
ably heavier computation than undertaken in this current 
work. 

In Fig. ((9]) we plot the lowest excitation frequencies 
as a function of temperature for the pancake trap with 
aspect ratio 1:7, (containing the above mentioned bi- 
concave structure). The m = 1 is as usual the Kohn 
mode. Note that at T = the m — 2 and m = 3 modes 
are lower than the m = 1 mode. The near degeneracy of 
the m = 1 and m = 2 modes at T = is accidental: for 
higher number of particles, the m = 2 and m — 3 modes 
actually goes below the m = 1 mode, a consequence of 
the discrete roton-like spectrum discussed in Ref. 
Close to the critical temperature the excitation energies 
approach their ideal gas values, as might be expected due 
to the decreasing density of the gas. In between there is 
an interesting crossing between the m = and m = 3 
modes. 
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D. Comparison with Monte Carlo 

Finally, we have attempted to compare the HFB-Popov 
method with the path integral Monte Carlo simulations 
of Nho and Landau @ . However we do not find a good 
agreement. The energies with HFB-Popov approxima- 
tion for a dipolar condensate in a pancake trap came 
about 10% higher than the Monte Carlo simulation, with 
a similar discrepancy in the shape (width) of the dipolar 
condensates. We note that the Monte Carlo simulations 
were performed for a very small number of particles, be- 
tween 27 and 125. Under these conditions, the critical 
temperature is very small, of the order of the trap fre- 
quencies, and in fact lower than the chemical potential. 
Thus, only a few low modes are excited even close to 
the critical temperature T c . At the simulated temper- 
ature of 0.4T C the density of the thermal gas is of the 
order of that of the condensate and there is a large over- 
lap between the two. Under these conditions it may be 
expected that the approximation of ignoring the thermal- 
thermal dipolar exchange interactions is invalid. Thus, it 
is plausible that the disagreement is due to this additional 
approximation rather than the inadequacy of the HFB- 
Popov method. But, for a number of particles of order 
10 4 — 10 5 , the critical temperature is much higher than 
the trap frequencies and the density of the thermal cloud 
significantly lower than that of the condensate. Thus, we 
believe our method should still give valid results under 
normal experimental conditions. 

IV. CONCLUSIONS 

In conclusion, we applied the Hartree-Fock- 
Bogoliubov-Popov approximation to dipolar gases 



in harmonic, cylindrically symmetric traps. For com- 
putational reasons, the exchange interaction due to the 
thermal gas has to be ignored (i.e, we ignore exchange 
interaction due to long range spatial correlation in the 
thermal component). For normal configurations where 
the condensate structure at T = has a maximum at the 
center, we observe a temperature dependent behavior 
similar to that of a gases with contact interaction 
[IH Ojj], but the behavior depends on the aspect ratio 
of the trap: for pancake traps, the dipolar interaction is 
effectively repulsive, leading to reduction of condensed 
part at a given temperature, and thus also reduction of 
the critical temperature for condensation, while in cigar 
trap, it is effectively attractive, leading to the opposite 
thermal effects. For configurations where bi-concave 
structure exists at T = 0, we find the somewhat surpris- 
ing result that this structure is actually enhanced at low 
temperatures (i.e, the ratio of the central density to the 
maximal density is reduced). For higher temperatures, 
the bi-concave structure gradually becomes less distinct 
and it disappears at a temperature of about 75% of the 
critical temperature. The low excitation spectrum of the 
bi-concave structure also shows interesting temperature 
dependence with crossing between different modes. 
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